Universal properties of a trapped two-component Fermi gas at unitarity 
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We treat the trapped two-component Fermi system, in which unlike fermions interact through 
a two-body short-range potential having no bound state but an infinite scattering length. By 
accurately solving the Schrodinger equation for up to TV = 6 fermions, we show that no many- 
body bound states exist other than those bound by the trapping potential, and we demonstrate 
unique universal properties of the system: Certain excitation frequencies are separated by 2hu), the 
wavefunctions agree with analytical predictions and a virial theorem is fulfilled. Further calculations 
up to N — 30 determine the excitation gap, an experimentally accessible universal quantity, and it 
agrees with recent predictions based on a density functional approach. 

PACS numbers: 



Ultracold Fermi gases are pure and controllable sys- 
tems where few- and many-body properties can be 
studied. Experiments now routinely convert a weakly- 
interacting atomic Fermi gas into a molecular Bose- 
Einstein condensate by tuning the scattering length near 
a Fano-Feshbach resonance. While the limiting behav- 
iors of this BCS-BEC crossover are well understood, 
the strongly-interacting regime remains fundamental and 
challenging. In particular, the unitary regime, character- 
ized by a diverging two-body s-wave scattering length a s , 
has been the subject of much experimental and theoret- 
ical research in the cold-atom community (see Ref. [l| 
for a review and an exhaustive list of Refs.). Of course, 
two-component Fermi gases at unitarity also greatly in- 
terest the nuclear physics community. Cold atomic Fermi 
gases and neutron matter with large a s are two physi- 
cal realizations of the idealized system described in the 
"Bertsch problem" [3] • Despite its importance, few essen- 
tially exact calculations exist for its energy eigenstates, in 
part because the system is strongly-correlated and has no 
small parameter for controlled perturbative treatments. 

Here, we solve the many-body Schrodinger equation 
for the trapped unitary Fermi gas with an even or odd 
number of fermions by two different numerical tech- 
niques: a basis set expansion technique that uses corre- 
lated Gaussians (CG) and a fixed-node diffusion Monte 
Carlo (FN-DMC) method. Our energies provide much 
needed benchmarks for few-body systems with short- 
range interactions that are important for the atomic, nu- 
clear, and condensed matter physics communities. A key 
result of our study is a demonstration that this system 
exhibits unique universal properties resembling those of 
the non-interacting system. 

For the class of two-body potentials that support no s- 
wave bound state, we find that many-body bound states 
do not exist for systems with N < 6 fermions. Our cal- 
culated excited-state spectrum for N < 6 explicitly con- 
firms, to within our numerical accuracy, the assertion 
by Werner and Castin l _3| that the excitation spectrum 
consists of families of breathing-mode states separated 



by Ihijj. These excitations are non-trivial in the sense 
that they are not associated with center-of-mass (CM) 
excitations. We interpret this spacing using hyperspher- 
ical coordinates, and show that in addition to these fam- 
ilies of breathing-mode excited states, other excitation 
frequencies exist that equal non-integer multiples of the 
trapping frequency. Moreover, the wavefunctions agree 
with analytically determined solutions 0], whose func- 
tional form was derived on the basis of a universality 
assumption. Finally, we show that a virial theorem is 
fulfilled 0,0]. This, together with our finding that neg- 
ative energy states do not exist, shows explicitly that 
the entire spectrum of systems with N < 6 is univer- 
sal; previously this had only been shown explicitly for 
the three-body system [y] but its likely validity had 
also been speculated for larger N (see, e.g., Ref. [3|, and 
Refs. therein). A final key point of the present study is 
a calculation of the excitation gap, for TV < 29. Encour- 
aging agreement is found with results obtained within a 
density functional theory (DFT) framework [7j, and it 
also sheds light on the applicability of the local density 
approximation (LDA). 

Our starting point is the Hamiltonian H for N mass- 
m fermions under spherical external confinement with 
angular trapping frequency u>, H = T + V tr + Vi n t, 
where T = -ft 2 /(2m) £ti V? , V tr = Eti m ^/ 2 



and V int = Y,?=iY,k=iMrjk). Here /Vj and N 2 de- 
note respectively the number of spin-up and spin-down 
fermions (N = N1+N2), and r*; denotes the position vec- 
tor of the ith atom. The purely attractive short-range 
potentials Vq (identical to those of Ref. [^) depend on 
the interparticle distance r^fc, and are parametrized by 
two parameters, the depth d and the range Rq, where 
i?o <C aho and a^ D = ^/hj (mto). For a given i?o, d is ad- 
justed so that the first free-space two-body s-wave bound 
state is just about to exist (i.e., so that |et s | = 00). To 
study even-odd oscillations, we set iVi = N% for even N, 
and Ni — N2 + 1 for odd N. Two numerical methods 
are applied to solve the Schrodinger equation for H: the 
CG and the FN-DMC approach. 
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TABLE I: CG energies E vn for a Gaussian potential; Ro = 
O.Olaho for N = 3 and 4, and 0.05a^ o for N = 5 and 6. 
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The CG approach has previously been used to deter- 
mine the energy spectrum of the four-fermion system 
with short-range interactions 0, Q . Here we extend these 
studies up to N — 6, and additionally present structural 
properties. In our implementation, each basis function is 
written as a product of the CM ground state and a sym- 
metrized product of Gaussian functions each of which de- 
pends on one of the N(N — l)/2 interparticle distances. 
The resulting states have vanishing relative orbital an- 
gular momentum L re i. To describe the iV-atom states 
with arbitrary L re i in our CG approach, we add an extra 
non-interacting particle. The Hamiltonian matrix is eval- 
uated analytically. A stochastic variational approach [10| 
optimizes the basis set and convergence. 

Table |l] summarizes a few selected total energies E vn 
for N = 3 — 6 (throughout this work, CM excitations 
are not considered, i.e., Ecu = |^)- For N = 3 and 
4, calculations for different Ro indicate that the finite 
range effect of the reported energies is w 0.02/zw 0, Q. 
For N = 5 and 6, finite range effects are expected to be 
f=a 0.05fiw. Assignment of the quantum numbers v and n 
is discussed below. Our N = 3 energies agree with those 
reported in Ref. [f| and the corresponding s v coefficients 



(defined below) agree with those reported in Ref. 
For N = 5 and 6, our energies are the first ab initio 
results obtained for this Hamiltonian by a method that 
is at least in principle free of any assumptions. We find 
that the systems with N = 5 (L re i = and 1) and N = 6 
(L re i — 0) support no negative energy states. 

To treat up to N = 30 fermions, the Schrodinger equa- 
tion is solved by the FN-DMC method [H, El- The 
proper fermionic antisymmetry is imposed through the 
use of a guiding function ipT- To within statistical un- 
certainties, the FN-DMC algorithm provides an upper 
bound to the exact energy, i.e., to the lowest-lying state 
with the same symmetry as ifa. Following Ref. [8|, we 
consider two different functional forms for ipx '■ The nodal 
surface of ipTi is constructed by antisymmetrizing a prod- 
uct of pair functions 1J| , and that of ipT2 coincides with 
the nodal surface of the non-interacting Fermi gas. Odd 
N systems are studied using generalizations of the tpr 
defined in Ref. Q for even N: We find that -0T2 gives 
the lowest energy up to N < 9 and ipri for larger N. 



TABLE II: FN-DMC energies Eoo and expectation values 
2{Vtr), both in units of hw. The range Ro of the square-well 
potential is O.Olaho for TV < 20, and a bit smaller for larger N. 
Statistical uncertainties of Eoo are reported in round brackets 
and of 2{Vtr) are a few times larger than those of _Eoo- 
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Unlike the energy, the unbiased calculation by the FN- 
DMC approach of observables associated with operators 
that do not commute with the Hamiltonian is more in- 
volved. Here, we use the mixed estimator {A) m i xe d, 
{A) mixed = 2(A) D mc ~ (A) V mc [3 • To impro ve ip T , 
we introduce additional two-body correlations [3, [l6| . 

Table ITT1 summarizes the ground state energies -Eoo ob- 
tained by the FN-DMC approach for N < 30 (the en- 
ergies for even N, N < 20, were already presented in 
Ref. Q). For N < 22, our energies provide more strin- 
gent bounds with smaller error bars than those reported 
previously [17] . For N = 23 — 30, these are the first ab 
initio results available. Comparison with Table |l] shows 
that the FN-DMC energies for N < 6 agree with the 
CG energies to within 2%. This agreement validates our 
construction of the nodal surface of ipT for N < 6. 

Dilute gases with short-range interactions exhibit 
unique universal properties at unitarity, such as the pre- 
dicted behavior that the energy spectrum contains se- 
quences of excited universal breathing-mode states sep- 
arated by 2Huj d, [lij]. This unique feature, specific to 
the unitary gas, can be intuitively understood by real- 
izing that the only length scale of the problem is given 
in this limit by the size of the system, i.e., the hyperra- 
dius R. Here, R is defined by removing the CM vec- 
tor Rcm and dividing the remaining 3iV — 3 coordi- 
nates into the hyperradius R and 3A^ — 4 hyperangles tt: 
mR 2 = J2i mr i ~~ MRjq M , where M — Nra. A dimen- 
sionality argument then suggests that the functional form 
of the universal hyperradial potentials V Sv (R) should be 
the same as that of the non-interacting system, which 
immediately implies the 2hw energy spacing for any ex- 
citation operator that depends on only R. 

This argument has been formalized by Werner and 
Castin for the N particle unitary gas interacting through 
zero-range pseudopotentials under harmonic confine- 
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ment It has been predicted that the adiabatic ap- 
proximation is exact for universal states. For these 
states, the wavefunction thus separates into a product 
of an in general unknown i?-independent channel func- 
tion $„(f2) and a ^-independent radial function F un (R), 
*"'(#, ft) = R^- 3N ^ 2 F un (R)<P u (Q). For a given hy- 
pcrangular quantum number v (v = 0, 1, • • • ), the radial 
quantum number n takes the values n — 0, 1, • • • . The 
description of the strongly-interacting many-body system 
thus reduces to solving a one-dimensional Schrodingcr 
equation for the hyperradial potential V Sv (i?) , which in- 
cludes part of the kinetic energy and a contribution due 
to the two-body interactions, 



-h 2 d 2 

^dR 2 



mLU 2 R 2 



F rel F 



(1) 



Here, V Su (R) = h 2 s v {s v + \)/{2mR 2 ). The ener- 
gies E r v ^ are related to the total energies E vn through 
E un = E%jl + Ecm- The coefficients s v are con- 
stants that arise from the integration over the hyper- 
angular Schrodinger equation; their values are, in gen- 
eral, unknown. Owing to the simple functional depen- 
dence of Eq. {T]) on R, the E^ and F vn can be writ- 
ten down readily @|, E^ = (a v + 2n + 3/2) frw and 



F m (R) = R s » +1 L { n 



. + 1/2) 



(R 2 /a 



ho) ex P( 



-R 2 /2al ) (not 



normalized) , where n denotes a non-negative integer, and 
£(s„+i/2) j g L a g Uerre polynomial. The expression for 
the E^ reveals immediately that the spacing between 
energy levels within a given hyperradial potential curve 
V Sv {R) is independent of s u and equals 2huj. 

Using the expression for E^, the Eqq reported in Ta- 
bles U and [ill readily give so for N < 30. For N — ► oo, the 
ground state energy of the trapped and homogenous sys- 
tems can be connected via the LDA (see, e.g., Ref. [191), 
leading to sq « y 'ihomE 'n i / 'fu^ '■ Here, E^j denotes 
the energy of the non-interacting trapped system and 
(,hom a universal parameter of the bulk system, £hom = 
0.42 [3 EE H ED], which leads to s « 0.65E NI /fuv. 
For the trapped system, we obtained £ tr = 0.465 [8j and 
so ~ 0.68Eni /hiu. The functional form of V Sv can also be 
obtained using a renormalization technique as N — > oo 
[13], which leads to so « 0.7lE^j /Hlo in fairly close 
agreement with the above values. 

To date, the 2hu> spacing of the unitary gas has only 
been verified explicitly for N = 3 The spacing of 
the energies E vn reported in Table [I] for N = 3 — 6 for 
the first few states with v = agrees with the predicted 
2hw spacing to better than 2%, i.e., within our numerical 
uncertainty. We additionally verified for N < 4 that this 
spacing holds for the lowest states with v = 1 and 2. 

Dash-dotted lines in Figs. [2(a) and (b) show the 
potential curves V(R), V(R) = V Su {R) +muj 2 R 2 /2, for 
N = 4 for v = and v = 1, respectively, calculated 
using the s v coefficients obtained from the CG energies. 
Dash-dotted lines in Fig. [TJc) show V(R) for v = for 
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FIG. 1: (Color Online) Dash-dotted lines show V(R) as a 
function of R/a ho for (a) N = 4 with v — 0, (b) N = 4 with 
v = 1, and (c) N = 9, 12 and 17 (from bottom to top) with 
v = 0. Circles and solid lines show the corresponding F„ n de- 
termined numerically and analytically, respectively. Dashed 
horizontal lines show the energies E r J^. 



N — 9, 12 and 17 (from bottom to top) calculated using 
the so coefficients obtained from the FN-DMC energies. 
Circles and solid lines show the corresponding square of 
the radial functions F 2 n obtained numerically (by inte- 
grating ^ 2 n over all coordinates but the hyperradius R) 
and analytically (using the s u coefficients calculated from 
the CG and FN-DMC energies), respectively. Clearly, 
our numerically determined F 2 n agree well with those 
determined analytically, providing the first quantitative 
numerical verification of the interpretation of the 2ftw 
spacing within the hyperspherical framework for N > 3. 

Another property of the unitary Fermi gas is that 
all of its universal states follow a virial theorem, i.e., 




FIG. 2: (Color Online) Circles show the excitation gap A(JV) 
determined from our FN-DMC energies. A solid line shows 
Alda(N) using £ tr and a tr . For comparison, triangles show 
A(iV) determined from the DFT energies Q. 

E vn = 2{V t r)^i n [!, H- Table |TT] shows that the en- 
ergy expectation value (columns 2 and 5) and the energy 
obtained using the virial theorem (columns 3 and 6) for 
N < 20 agree well, providing further evidence for the 
separability of universal states of the unitary gas and, 
on the other hand, suggesting that our parametrization 
of the nodal surface used in the FN-DMC calculations is 
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adequate. The virial theorem has previously been tested 
at various levels of approximation 0, [U Hd] . 

We now use the ground state energies reported in Ta- 
ble [IT] to determine the experimentally measurable exci- 
tation gap A(iV) for the trapped unitary Fermi g clS clS ct 
function of N. The definition of the excitation gap for the 
homogeneous system [2(J, which equals half the energy 
it takes to break a pair, can be readily extended to the 
trapped system A (AT) = E 00 (N) - [E oa (N - 1) + 

Eoo(N + l)]/2. Circles in Fig. [2]show the excitation gap 
A(N) obtained from our FN-DMC energies. A(N) in- 
creases from w 0.75hiu for N = 3 to w 1.5ftw for A" = 29. 
For comparison, triangles show A(N) calculated using 
the DFT energies obtained recently by Bulgac Q- The 
agreement between our and Bulgac's A(iV) is quite good. 

To gain further insight, we calculate the gap in the 
LDA, A LDA {N) = 3a /lom (3A^) 1 /3^/(8V&). The 
universal parameter a describes the even-odd oscilla- 
tions: ahom equals 0.85 for the bulk system [2l[, and 
we find atr — 0.60 for the trapped system. Alda(N), 
shown in Fig. [2] by a solid line using £tr and atr, provides 
a good description of our FN-DMC results. The fact 
that atr is noticeably smaller than ahom suggests that 
the extra particle is not distributed uniformly through- 
out the cloud. Indeed, our density profiles for TV > 11 
(not shown) indicate that the extra particle sits near the 
surface of the cloud. In the LDA, the center should be 
described well but not the surface, so LDA might fail 
here. Recently, Son proposed that the gap increases with 



24j. A fit of our results for A(N) for 



N 1 / 9 as N 

N > 9 shows consistency with the N 1 / 9 dependence but 
does not conclusively confirm it. 

In summary, we have shown that trapped two- 
component unitary Fermi gases with N < 6 support no 
negative energy states. While our FN-DMC calculations 
do not exclude the presence of negative energy states for 
larger N, we find no evidence for their existence. A vari- 
ational analysis predicts an upper bound for the system 
size N cr at which non-universal states exist of « 40. To 
obtain N cr we consider the Gaussian two-body potential 
and take the variational many-body wave function of the 
untrapped system to coincide with that of a trapped non- 
interacting two-component Fermi gas. The confinement 
width is treated as a variational parameter, and N is in- 
creased till a negative energy state is found. The absence 
of negative energy states for N < 6 fermions explains the 
stability of equal-mass two-component Fermi gases. Fur- 
thermore, we explicitly verify a number of unique prop- 
erties of the universal states of the unitary Fermi gas. 
The particular form of the energy spectrum has conse- 
quences: It should be possible to verify experimentally 
that for each excitation frequency uj ex there exists a lad- 
der of non CM excitations uj ex +2nuj (n integer). Further- 
more, the s v coefficients calculated here for the trapped 



system remain valid for the free-space system Q. Fi- 
nally, we calculate the excitation gap A(N) for A" < 29, 
and find good agreement with recent DFT results. A (AT) 
increases with A?" and the excess atom leads for A^ > 11 
to an increase of the density in the surface region. 
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